.

mae <- readRDS("~/BHK lab/Ravi/Ravi_Test/data/result3.rds")

#removing na in response and expression data

na_rows <- is.na(mae$response)
cleaned_response <- mae$response[!na_rows]
cleaned_expr_data <- assays(mae)[["expr"]][, !na_rows]


# 122 patient ate removing NA
design <- model.matrix(~ cleaned_response)
fit <- lmFit(cleaned_expr_data, design)
fit <- eBayes(fit)

#Top table only shows the top portion of the results
topTable(fit)
## Removing intercept from test coefficients
##                       logFC   AveExpr        t      P.Value adj.P.Val
## ENSG00000211835.1 2.9729828 -7.868755 4.290500 3.620134e-05 0.9693367
## ENSG00000173231.6 1.5775015 -8.827553 3.826189 2.077302e-04 0.9693367
## ENSG00000229037.2 1.5434391 -8.867183 3.734687 2.884608e-04 0.9693367
## ENSG00000261335.1 1.0303203 -9.245205 3.693547 3.337526e-04 0.9693367
## ENSG00000251497.2 1.0942021 -9.070539 3.660129 3.754202e-04 0.9693367
## ENSG00000244053.1 1.7543381 -2.320355 3.491816 6.713542e-04 0.9693367
## ENSG00000236651.1 0.8013953 -9.597931 3.478346 7.027426e-04 0.9693367
## ENSG00000232549.1 2.0167536 -7.511261 3.469648 7.237379e-04 0.9693367
## ENSG00000232173.2 1.4258124 -8.763200 3.465522 7.339043e-04 0.9693367
## ENSG00000266707.1 2.2187534 -6.882755 3.457535 7.539641e-04 0.9693367
##                             B
## ENSG00000211835.1  1.85507430
## ENSG00000173231.6  0.41670113
## ENSG00000229037.2  0.14769222
## ENSG00000261335.1  0.02836268
## ENSG00000251497.2 -0.06781811
## ENSG00000244053.1 -0.54187398
## ENSG00000236651.1 -0.57905328
## ENSG00000232549.1 -0.60299910
## ENSG00000232173.2 -0.61434367
## ENSG00000266707.1 -0.63627069
datatable(topTable(fit))
## Removing intercept from test coefficients

Preparing Data for volcano Plot

Convert fit to Data frame and add column gene_name symbol

result <- topTable(fit, number=Inf)  # Get all results
## Removing intercept from test coefficients
df <- as.data.frame(result)  # Convert to data frame

#extract gene_name from mae ( it has gene_id column and gene_name)

genedata<- data.frame(rowData(mae@ExperimentList$expr))

#subset gene_id and gene_name from the gene_data

subset_genedata <- genedata[, c("gene_id","gene_name" )]

#mergeing process

#add a gene_id column to result 
result$gene_id <- rownames(result)
#merge result and subset_genedata  by gene_id 

merge_result <- merge(result, subset_genedata, by= "gene_id")

datatable(merge_result)
## Warning in instance$preRenderHook(instance): It seems your data is too big for
## client-side DataTables. You may consider server-side processing:
## https://rstudio.github.io/DT/server.html

Volcano Plot

EnhancedVolcano(merge_result,
    lab = merge_result$gene_name,
    x = 'logFC',
    y = 'P.Value')